Aminet 1 (Walnut Creek)/Aminet - June 1993 [Walnut Creek].iso/usenet/sources/volume91/utilitys/ed3_1_1/part06/edvan.c

C/C++ Source or Header | 1991-11-07 | 32KB | 1,474 lines

/*
    These are miscellaneous vanilla functions of ed3,
        using no graphics or intuition calls.

10-22-91 chg YZ to ZY
         fixed a very small (but serious) bug in find_dual
10-13-91 added fibo() fibonacci (?) function
         stellate() now checks to see if there is room in the undo buf
09-25-91 chg points_are_a_face2 to face_next_to_edge
         add adjacent_fp(), adjacent_pp(), adjacent_xx()
         fix a problem in find_dual with polys->polys
09-06-91 adapted to use coord instead of short
         rem distance_between_points
09-04-91 chg face_next_to_edge() to allow an ignored face
09-03-91 chg rotate() to support the rotate_all flag
         implement deselect_connected()
09-02-91 add center_to_origin(), move_point_rel(), move_point_to()
08-27-91 add points_are_on_a_poly()
08-26-91 fix to point_is_on_a_face()
         add point_is_on_a_poly()
08-24-91 add stellate() (hey, it was easy!)
         add glue_points()    (also easy)
08-23-91 added refresh_all(), toggle_dtool()
08-21-91 added select_all(), select_connected(), deselect_all(),
         deselect_connected(), toggle_select();
08-15-91 added find_center()
         added find_radius()
         chg find_dual to support assume_poly flag
*/

#include "sysnogr.h"


#if INTEGER
void center_to_origin()
{
    double dcx, dcy, dcz;
    pixel cx, cy, cz;
    index p;

    if (!points) return;

    find_center(&dcx, &dcy, &dcz);
    cx = dcx;
    cy = dcy;
    cz = dcz;

    if (!cx && !cy && !cz)
    {
        screen_say("Points are already centered");
        return;    /* already at center */
    }

    for (p = 0; p < points; p++)
      move_point_rel(p, -cx, -cy, -cz);

    refresh_all();
}
#else
void center_to_origin()
{
    double dcx, dcy, dcz;
    index p;

    if (!points) return;

    find_center(&dcx, &dcy, &dcz);

    if ((dcx > -.001 && dcx < .001) &&
        (dcy > -.001 && dcy < .001) &&
        (dcz > -.001 && dcz < .001))
    {
        screen_say("Points are already centered");
        return;    /* already at center */
    }

    for (p = 0; p < points; p++)
      move_point_rel(p, -dcx, -dcy, -dcz);

    refresh_all();
}
#endif


void reset_zoom()
{
  SCALE = 1.0;
  offset.x = offset.y = offset.z = 0;
}


void move_point_rel(p, dx, dy, dz)
index p;
coord dx, dy, dz;
{
    action_move_p(p, dx, dy, dz);

    point[p].x += dx;
    point[p].y += dy;
    point[p].z += dz;
    convert_point(p);
}


void move_point_to(p, x, y, z)
index p;
coord x, y, z;
{
    action_move_p(p, x-point[p].x, y-point[p].y, z-point[p].z);

    point[p].x = x;
    point[p].y = y;
    point[p].z = z;
    convert_point(p);
}


void refresh_all()
{
    int i1;

    draw_faces(1);
    draw_points(0);
    if (flags.display_dtool)
        draw_distance_tool();
    if (flags.show_coords) draw_coords(0);
    if (mode == M_MORE || mode == M_KMORE)
        redraw_newface();    /* redraw the partially completed face */
    if (mode == M_MORE)
    {
        if (flags.shape == M_FACE)
            for (i1 = 0; i1 < new_seg; i1++) tri_up_point(face[faces].p[i1]);
        else
            for (i1 = 0; i1 < new_seg; i1++) tri_up_point(temp_poly[i1]);
    }
    if (mode == M_KMORE)
    {
        if (flags.shape == M_FACE)
            for (i1 = 0; i1 < new_seg; i1++) tri_down_point(face[faces].p[i1]);
        else
            for (i1 = 0; i1 < new_seg; i1++) tri_down_point(temp_poly[i1]);
    }
    if (flags.crosshair) draw_hair();            /* Erase/Redraw crosshair */
    draw_bars();
    show_counts();
}


void toggle_dtool()
{
    flags.display_dtool ^= 1;
    draw_distance_tool();
    /* if we turn the dtool off, exit dtool mode */
    if (!flags.display_dtool && mode == M_DIST)
    {
        mode = M_ADDP;
        draw_bars();
    }
}


void select_all()
{
    index p, limit;

    deselect_all();

    limit = (points >= MAX_SELECT ? MAX_SELECT : points);
    for (p = 0; p < limit; p++)
    {
        select_p[p] = p;
        box_point(p);                /* select it visually */
        point[p].code |= PF_SELECTED;    /* set selected flag */
    }
    selected = p;
}


void deselect_all()
{
    index i1;

    for (i1 = 0; i1 < selected; i1++)
    {
        point[select_p[i1]].code &= ~PF_SELECTED;    /* turn off selected flag */
        box_point(select_p[i1]);                    /* and deselect it visually */
    }
    selected = 0;
}


/*
    If the starting point is specified, the search limits itself to points
    connected to the starting point.  Otherwise, the search starts from
    every selected point
*/

void select_connected(starting_point)
index starting_point;
{
    int sp, i1, i2, i3, verts;    /* selected point, face and polygon indices */
    index start_p, p;
    index orig_selected = selected;

    start_p = 0;

    if (starting_point >= 0)
    {
        for (i1 = 0; i1 < selected; i1++)
            if (starting_point == select_p[i1])
            {
                start_p = i1;
                break;
            }
    }

    for (sp = start_p; sp < selected; sp++)
    {
        p = select_p[sp];

        for (i1 = 0; i1 < faces; i1++)
        for (i2 = 0; i2 < 3; i2++)
            if (face[i1].p[i2] == p)
            {
                /* select the other points of this face */
                select_point(face[i1].p[(i2+1)%3]);
                select_point(face[i1].p[(i2+2)%3]);
            }

        for (i1 = 0; i1 < polys; i1++)
        {
            verts = poly[i1].verts;
            for (i2 = 0; i2 < verts; i2++)
            {
                if (poly[i1].p[i2] == p)
                {
                    /* select the other points of this poly */
                    for (i3 = 1; i3 < verts; i3++)
                    select_point(poly[i1].p[(i2 + i3) % verts]);
                }
            }
        }

        if (starting_point >= 0 && sp == start_p)
            sp = orig_selected-1;
    }
}


#define DS_MAX 500
index dselect[DS_MAX];
index dselected;

void deselect_connected(starting_point)
index starting_point;
{
    int i1, i2, i3, verts;    /* selected point, face and polygon indices */
    index dp, p, newp;

    dselected = 1;
    dselect[0] = starting_point;

    for (dp = 0; dp < dselected; dp++)
    {
        p = dselect[dp];

        for (i1 = 0; i1 < faces; i1++)
        for (i2 = 0; i2 < 3; i2++)
        if (face[i1].p[i2] == p)
        {
            /* deselect the other points of this face */
            newp = face[i1].p[(i2+1)%3];
            if (pselected(newp))
            {
                deselect_point(newp);
                dselect[dselected++] = newp;
                if (dselected == DS_MAX) return;
            }

            newp = face[i1].p[(i2+2)%3];
            if (pselected(newp))
            {
                deselect_point(newp);
                dselect[dselected++] = newp;
                if (dselected == DS_MAX) return;
            }
        }

        for (i1 = 0; i1 < polys; i1++)
        {
            verts = poly[i1].verts;
            for (i2 = 0; i2 < verts; i2++)
            if (poly[i1].p[i2] == p)
            for (i3 = 1; i3 < verts; i3++)    /* select the other points of this poly */
            {
                newp = poly[i1].p[(i2 + i3) % verts];
                if (pselected(newp))
                {
                    deselect_point(newp);
                    dselect[dselected++] = newp;
                    if (dselected == DS_MAX) return;
                }
            }
        }
    }
}


void toggle_select(buttons, tp)
int buttons;
index tp;
{
    if (buttons & 2)    /* double click: select connected */
    {
        if (pselected(tp))
            select_connected(tp);
        else
            deselect_connected(tp);
    }
    else        /* toggle individual point */
    {
        if (pselected(tp))
            deselect_point(tp);
        else
            select_point(tp);
    }
}


void glue(nearest)
index nearest;
{
    index second;

    if (nearest < 0) return;    /* must have valid point */
    second = find_nearest_in_space(point[nearest].x, point[nearest].y,
            point[nearest].z, nearest);
    if (second < 0) return;
    begin_operation();
    glue_points(nearest, second);
    end_operation();
    refresh_all();
}


/*
    glue_points fuses two points together into one point.
    They should be spatially close to begin with.
*/

void glue_points(p0, p1)
{
    Point3D new;    /* coordinate of new fused point */
    index i1, vert;
    index changed[MAX_EDGES];
    index found = 0;
    index f, q;
    index fp0, fp1, fp2;

    new.x = (point[p0].x + point[p1].x)/2;
    new.y = (point[p0].y + point[p1].y)/2;
    new.z = (point[p0].z + point[p1].z)/2;

    /* handle possible degenerations of faces and polys */
    for (f = 0; f < faces; f++)        /* search all faces */
    {
        fp0 = face[f].p[0];
        fp1 = face[f].p[1];
        fp2 = face[f].p[2];

        if ((fp0 == p0 && fp1 == p1) ||
            (fp1 == p0 && fp0 == p1) ||
            (fp0 == p0 && fp2 == p1) ||
            (fp2 == p0 && fp0 == p1) ||
            (fp1 == p0 && fp2 == p1) ||
            (fp2 == p0 && fp1 == p1))
        {
            del_face(f, 0);    /* faces cannot degenerate into lines */
            f--;    /* don't skip a face */
        }
    }
    /* now look for a polygon having these two points as an edge */
    /* yeah, calling a function here makes things slower, but hey, speed is
        probably not important.  If it was, we could use our own loop here. */
    while ((q = poly_next_to_edge(p0, p1, -1)) != -1)
    {
        /* which vertex is it? */
        for (vert = 0; vert < poly[q].verts; vert++)
            if (poly[q].p[vert] == p1)
            {
                delete_vertex(q, vert, 0);
                break;    /* breaks out of for loop */
            }
        if (poly[q].verts == 3)    /* we have reduced it to a triangle */
        {
            poly_to_face(q, 0);
            del_poly(q, 0);
        }
    }

    /* change all instances of the second point into the first */
    for (i1 = 0; i1 < faces; i1++)
        for (vert = 0; vert < 3; vert++)
            if (face[i1].p[vert] == p1)
            {
                changed[found++] = i1;
                face[i1].p[vert] = p0;
            }

    for (i1 = 0; i1 < polys; i1++)
        for (vert = 0; vert < poly[i1].verts; vert++)
            if (poly[i1].p[vert] == p1)
            {
                changed[found++] = i1 | POLYFLAG;
                poly[i1].p[vert] = p0;
            }
    /* the join_p event is really only necessary if there were facets
        referencing the eliminated point */
    if (found) action_join_p(p0, p1, changed, found);

    del_point(p1, 1, 0);    /* kill quickly: assume the point is loose */
                            /* but register the action */
    move_point_to(p0, new.x, new.y, new.z);
}


void unfold()
{
}


void stellate()
{
    int go = 1, event;
    int id;
    long data = 100, mag = 100;    /* magnitude of stellation */
    int rel_to;
    int add_p, add_f, del_f, del_q;
    index po0;
    long total_verts = 0;    /* total edges of all polygons */
    int max_verts = 0;        /* the most verts on any polygon */
    int verts;

    sys_window(35, 4, "Stellate All Facets");
    sys_movecur(1,0);
    sys_say_text("Magnitude:");
    sys_get_long(data, 1);
    sys_movecur(1,1);
    sys_say_text("Stellate relative to the object's");
    sys_movecur(1,2);
    sys_boolean("Face Center", 2);
    sys_say_text(" or ");
    sys_boolean("Radius", 3);
    sys_activate(1);

    while (go)
    {
        event = sys_read_event(&id, &data);
        switch (event)
        {
        case EV_CANCEL:
            sys_close_window();
            return;
            break;
        case EV_CHECKINT:
        case EV_HITRETURN:
        case EV_HITBUTTON:
            if (id == 1)
            {
                if (data > 1000) data = 1000;
                if (data < -1000) data = -1000;
                mag = data;
            }
            else if (id == 2)
            {
                rel_to = 0;
                go = 0;
            }
            else if (id == 3)
            {
                rel_to = 1;
                go = 0;
            }
            break;
        }
    }
    sys_close_window();

    /* Estimate how many actions the stellation will require. */
    for (po0 = 0; po0 < polys; po0++)
    {
        verts = poly[po0].verts;
        total_verts += verts;
        if (max_verts < verts) max_verts = verts;
    }
    add_p = faces + polys;
    add_f = (faces * 3);        /* for the faces */
    add_f += total_verts;        /* for the polys */
    del_f = faces;
    del_q = polys;

    /* Make sure that this event won't overflow the undo buffer */
    if (!enough_room(add_p, add_f, 0, 0, del_f, del_q, 0, max_verts))
    {
        /* this operation will overflow the undo buffer! */
        if (!tell_user_overflow()) return;
    }
    begin_operation();
    stellate0(mag, rel_to);
    end_operation();
}


void stellate0(stellar_mag, rel_to)
int stellar_mag;    /* magnitude of the stellation */
int rel_to;            /* 0 = relative to face center, 1 = relative to radius */
{
    uindex
        p0, p1, i, j,
        old_points, old_faces, old_polys,
        edges;
    index np;            /* new point (the one we've just created) */
#if INTEGER
    long vx, vy, vz;
#else
    double vx, vy, vz;
#endif
    double cx, cy, cz, radius;    /* center coordnates for polyhedron */
    double new_radius, scale_factor;
    index new_faces;
    char text[80];

    if (!faces && !polys) return;    /* must have something to stellate */

    begin_operation();

    /* estimate how many new faces will be generated */
    new_faces = (faces * 3);
    for (i = 0; i < polys; i++)
        new_faces += poly[i].verts;

    if (faces + new_faces >= MAX_FACES)    /* need to make room for the new faces */
    {
        if (more_faces(faces + new_faces + 100))
        {
            sprintf(text, "increased face space to %d", MAX_FACES);
            screen_say(text);
        }
        else
        {
            sprintf(text, "not enough memory for %d new faces", new_faces);
            screen_say(text);
            return;                /* failure */
        }
    }
    old_points = points;    /* After we've created the stellated points */
    old_faces = faces;        /* and faces, we'll remove the old faces. */
    old_polys = polys;        /* And polygons. (but NOT the old points) */

    /* assume it's a polyhedron that we're operating on */
    find_center(&cx, &cy, &cz);
    find_radius(cx, cy, cz, &radius);

    /**** First we create the new points from the old faces ****/
    for (i = 0; i < old_faces; i++)
    {
        vx = vy = vz = 0;
        for (j = 0; j < 3; j++)
        {
            p0 = face[i].p[j];
            vx += (point[p0].x - cx);
            vy += (point[p0].y - cy);
            vz += (point[p0].z - cz);
        }
        /* project point at face center out to stellar length */
        vx /= 3;
        vy /= 3;
        vz /= 3;
        new_radius = sqrt((double)(vx*vx + vy*vy + vz*vz));
        if (rel_to)
            scale_factor = (new_radius + stellar_mag) / new_radius;
        else
            scale_factor = (radius + stellar_mag) / new_radius;
        np = add_point((coord)(cx + (double)vx * scale_factor),
                      (coord)(cy + (double)vy * scale_factor),
                      (coord)(cz + (double)vz * scale_factor));

        /* Add stellation */
        add_face(face[i].p[0], face[i].p[1], np, colors.newface, 0);
        add_face(face[i].p[1], face[i].p[2], np, colors.newface, 0);
        add_face(face[i].p[2], face[i].p[0], np, colors.newface, 0);
    }

    /* and from the old polygons */
    for (i = 0; i < old_polys; i++)
    {
        vx = vy = vz = 0;
        edges = poly[i].verts;
        for (j = 0; j < edges; j++)
        {
            p0 = poly[i].p[j];
            vx += (point[p0].x - cx);
            vy += (point[p0].y - cy);
            vz += (point[p0].z - cz);
        }
        vx /= edges;
        vy /= edges;
        vz /= edges;
        new_radius = sqrt((double)(vx*vx + vy*vy + vz*vz));
        if (rel_to)
            scale_factor = (new_radius + stellar_mag) / new_radius;
        else
            scale_factor = (radius + stellar_mag) / new_radius;
        np = add_point((coord)(cx + (double)vx * scale_factor),
                      (coord)(cy + (double)vy * scale_factor),
                      (coord)(cz + (double)vz * scale_factor));
        /* Add stellation */
        for (j = 0; j < edges; j++)
        {
            p0 = poly[i].p[j];
            p1 = poly[i].p[(j+1)%edges];
            add_face(p0, p1, np, colors.newface, 1);
        }
    }

    /*** Now we must destroy the original points, faces, and polys ***/
    for (i = 0; i < faces-old_faces; i++) face[i] = face[old_faces+i];

    /* make sure we release the memory held by the original polys */
    for (i = 0; i < old_polys; i++) kill_poly(i);
    for (i = 0; i < polys-old_polys; i++) poly[i] = poly[old_polys+i];

    faces -= old_faces;
    polys -= old_polys;

    end_operation();
}


void find_dual()
{
    uindex
        p0, i, j, k,
        af, adj_f[MAX_EDGES],    /* number and indices of the adjacent faces */
        ap, adj_p[MAX_EDGES],    /* number and indices of the adjacent polys */
        np, new_p[MAX_EDGES],
        old_points, old_faces, old_polys,
        edges;
    index used_f[MAX_EDGES], used_p[MAX_EDGES];
#if INTEGER
    long vx, vy, vz;
#else
    double vx, vy, vz;
#endif
    double cx, cy, cz, radius;    /* center coordnates for polyhedron */
    double new_radius, scale_factor;
    int flag;
#define DF_POLY 1
#define DF_FACE 2
    index previous;

    if (!faces && !polys) return;    /* must have something to dualize */
    if (!dual_ok()) return;

    begin_operation();

    old_points = points;    /* After we've created the dual's points */
    old_faces = faces;        /* and faces, we'll remove the old ones. */
    old_polys = polys;        /* And polygons. */

    if (flags.assume_poly)    /* assume it's a polyhedron that we're operating on */
    {
        find_center(&cx, &cy, &cz);
        find_radius(cx, cy, cz, &radius);
    }

    /**** First we create the new points from the old faces ****/
    for (i = 0; i < old_faces; i++)
    {
        vx = vy = vz = 0;
        for (j = 0; j < 3; j++)
        {
            p0 = face[i].p[j];
            if (flags.assume_poly)
            {
                vx += (point[p0].x - cx);
                vy += (point[p0].y - cy);
                vz += (point[p0].z - cz);
            }
            else
            {
                vx += point[p0].x;
                vy += point[p0].y;
                vz += point[p0].z;
            }
        }
        if (flags.assume_poly)    /* project point at face center out to radius length */
        {
            new_radius = sqrt((double)(vx*vx + vy*vy + vz*vz));
            scale_factor = radius / new_radius;
            add_point((coord)(cx + (double)vx * scale_factor),
                      (coord)(cy + (double)vy * scale_factor),
                      (coord)(cz + (double)vz * scale_factor));
        }
        else
            add_point(vx/3, vy/3, vz/3);    /* add a point at the face's center */
    }

    /* and from the old polygons */
    for (i = 0; i < old_polys; i++)
    {
        vx = vy = vz = 0;
        edges = poly[i].verts;
        for (j = 0; j < edges; j++)
        {
            p0 = poly[i].p[j];
            if (flags.assume_poly)
            {
                vx += (point[p0].x - cx);
                vy += (point[p0].y - cy);
                vz += (point[p0].z - cz);
            }
            else
            {
                vx += point[p0].x;
                vy += point[p0].y;
                vz += point[p0].z;
            }
        }
        if (flags.assume_poly)    /* project point at face center out to radius length */
        {
            new_radius = sqrt((double)(vx*vx + vy*vy + vz*vz));
            scale_factor = radius / new_radius;
            add_point((coord)(cx + (double)vx * scale_factor),
                      (coord)(cy + (double)vy * scale_factor),
                      (coord)(cz + (double)vz * scale_factor));
        }
        else add_point(vx/edges, vy/edges, vz/edges);    /* add a point at the poly's center */
    }

    /* then for each point, we find the surrounding faces */
    /* and make a face out of their corresponding points */

    for (i = 0; i < old_points; i++)
    {
        if (flags.debug) printf("Point %d: ", i);
        af = ap = 0;

        for (j = 0; j < old_faces; j++)
            if (face[j].p[0] == i || face[j].p[1] == i || face[j].p[2] == i)
                adj_f[af++] = j;

        for (j = 0; j < old_polys; j++)
            for (k = 0; k < poly[j].verts; k++)
                if (poly[j].p[k] == i) adj_p[ap++] = j;

        if (flags.debug)
        {
            if (flags.debug) printf("%d adjacent faces:", af);
            if (af) for (j = 0; j < af; j++) printf(" %d", adj_f[j]);
            printf("\n");
            printf("         %d adjacent polys:", ap);
            if (ap) for (j = 0; j < ap; j++) printf(" %d", adj_p[j]);
            printf("\n");
        }

        if (af+ap == 3)        /* three points are a simple face */
        {
            if (flags.debug) printf(" %d %d  adding face:", af, ap);
            np = 0;
            for (j = 0; j < af; j++)
                new_p[np++] = adj_f[j];
            for (j = 0; j < ap; j++)
                new_p[np++] = old_faces + adj_p[j];
            if (flags.debug)
            {
                printf("   ->");
                for (j = 0; j < 3; j++) printf(" %d", new_p[j]);
                printf("\n");
            }
            add_face(new_p[0], new_p[1], new_p[2], colors.newface, 1);
        }
        else if (af+ap > 3)
        {
            if (flags.debug) printf(" %d %d  adding POLYGON:", af, ap);
            for (j = 0; j < af; j++) used_f[j] = 0;
            for (j = 0; j < ap; j++) used_p[j] = 0;
            allocate_poly(af+ap, colors.newface);

            /*** At this point we must arrange the points in the correct
                order to produce a correct polygon. ***/
            /* We pick one arbitrary starting face, then advance around
                the point by finding adjacent faces/polys that we haven't
                used before. */
            np = 0;
            if (af)    /* pick starting point at 0 */
            {
                previous = poly[polys].p[np++] = adj_f[0];
                used_f[0] = 1;
                flag = DF_FACE;
            }
            else
            {
                previous = adj_p[0];
                poly[polys].p[np++] = old_faces + adj_p[0];
                used_p[0] = 1;
                flag = DF_POLY;
            }
            np = 1;
            while (np < af+ap)
            {
                for (k = 0; k < af; k++)    /* search faces */
                {
                    if (used_f[k]) continue;    /* skip already-used faces */
                    if (adjacent_xx(previous, adj_f[k], (flag == DF_POLY), 0))
                    {
                        previous = poly[polys].p[np++] = adj_f[k];
                        used_f[k] = 1;
                        flag = DF_FACE;
                    }
                }
                for (k = 0; k < ap; k++)    /* search polys */
                {
                    if (used_p[k]) continue;    /* skip already-used polys */
                    if (adjacent_xx(previous, adj_p[k], (flag == DF_POLY), 1))
                    {
                        poly[polys].p[np++] = old_faces + adj_p[k];
                        previous = adj_p[k];
                        used_p[k] = 1;
                        flag = DF_POLY;
                    }
                }
            }
            finish_poly(0);
            if (flags.debug)
            {
                printf("   ->");
                for (j = 0; j < af+ap; j++) printf(" %d", poly[polys-1].p[j]);
                printf("\n");
            }
        }
    }

    if (flags.debug)
    {
        printf("points %d, faces %d, polys %d\n", points, faces, polys);
        printf("starting destruction:\n");
    }

    /*** Now we must destroy the original points, faces, and polys ***/

    for (i = 0; i < points-old_points; i++)
        point[i] = point[old_points+i];
    for (i = 0; i < faces-old_faces; i++) face[i] = face[old_faces+i];

    /* make sure we release the memory held by the original polys */
    for (i = 0; i < old_polys; i++) kill_poly(i);
    for (i = 0; i < polys-old_polys; i++) poly[i] = poly[old_polys+i];

    points -= old_points;
    faces -= old_faces;
    polys -= old_polys;

    if (flags.debug)
    {
        printf("ending destruction:\n");
        printf("points %d, faces %d, polys %d\n", points, faces, polys);
    }
    convert_points();
    end_operation();
}


/* returns whether it is OK to perform a dualize operation */

int dual_ok()
{
    index p0, f0, q0, v0, surround;
    index p, f, q, dp, df, dq;

    if (!flags.undo_active) return 1;    /* always OK if undo is off */
    p = f = q = 0;

    /* Estimate how many faces and polys will need to be created */
    for (p0 = 0; p0 < points; p0++)
    {
        surround = 0;
        for (f0 = 0; f0 < faces; f0++)
            if (face[f0].p[0] == p0 ||face[f0].p[2] == p0 ||face[f0].p[1] == p0)
                surround ++;
        for (q0 = 0; q0 < polys; q0++)
            for (v0 = 0; v0 < poly[q0].verts; v0++)
                if (poly[q0].p[v0] == p0)
                    surround ++;
        if (surround == 3) f++;
        if (surround > 3) q++;
    }
    p = faces + polys;    /* each facet becomes a point */
    dp = points;
    df = faces;
    dq = polys;
    /* Make sure that this event won't overflow the undo buffer */
    if (!enough_room(p, f, q, dp, df, dq, 0, 0))
    {
        /* this operation will overflow the undo buffer! */
        if (!tell_user_overflow()) return 0;
    }
    return 1;
}


/*
    Compare a face to face, face to poly, poly to face, or poly to poly.
    The two flag arguments indicate what type of face the first two arguments
    are (flag = TRUE if they indicated argument is a poly)
*/

index adjacent_xx(qf0, qf1, poly_flag0, poly_flag1)
uindex qf0, qf1;
int poly_flag0, poly_flag1;
{
  if (poly_flag0)
  {
    if (poly_flag1)
        return adjacent_pp(qf0, qf1);
    else
        return adjacent_fp(qf1, qf0);
  }
  else
  {
    if (poly_flag1)
        return adjacent_fp(qf0, qf1);
    else
        return adjacent_ff(qf0, qf1);
  }
}


/*
    Compare two faces to see if they share an edge: returns 1 or -1
    depending on the direction of the match, 0 if there is no shared
    edge.
*/

index adjacent_ff(f0, f1)
uindex f0, f1;
{
    uindex e, e2, p00, p01, p10, p11;
    index *fpp = face[f1].p;

    for (e = 0; e < 3; e++)
    {
      p00 = face[f0].p[e];
      p01 = face[f0].p[(e+1)%3];
      for (e2 = 0; e2 < 3; e2++)
      {
        /* There are two ways that one edge can match another:
            the point numbers can be identical or inverted.
            Edge i has indices i, (i+1)%3 for i = 0,1,2.
        */
        p10 = fpp[e2];
        p11 = fpp[(e2+1)%3];

        if (p00 == p10 && p01 == p11) return 1;
        else if (p00 == p11 && p01 == p10) return -1;
      }
    }
    return (index)0;
}


/*
    Compare a face to a poly to see if they share any edge: returns 1 or -1
    depending on the direction of the match, 0 if there is no shared
    edge.
*/

index adjacent_fp(f0, q0)
uindex f0, q0;
{
    uindex e, e2, p00, p01, p10, p11, verts;
    index *fpp = face[f0].p;
    index *qpp = poly[q0].p;

    verts = poly[q0].verts;

    for (e = 0; e < 3; e++)
    {
      p00 = fpp[e];
      p01 = fpp[(e+1)%3];
      for (e2 = 0; e2 < verts; e2++)
      {
        p10 = qpp[e2];
        p11 = qpp[(e2+1)%verts];

        if (p00 == p10 && p01 == p11) return 1;
        else if (p00 == p11 && p01 == p10) return -1;
      }
    }
    return (index)0;
}


/*
    Compare a poly to a poly to see if they share any edge: returns 1 or -1
    depending on the direction of the match, 0 if there is no shared
    edge.
*/

index adjacent_pp(q0, q1)
uindex q0, q1;
{
    uindex e, e2, p00, p01, p10, p11, verts0, verts1;
    index *qp0 = poly[q0].p;
    index *qp1 = poly[q1].p;

    verts0 = poly[q0].verts;
    verts1 = poly[q1].verts;

    for (e = 0; e < verts0; e++)
    {
      p00 = qp0[e];
      p01 = qp0[(e+1)%verts0];
      for (e2 = 0; e2 < verts1; e2++)
      {
        p10 = qp1[e2];
        p11 = qp1[(e2+1)%verts1];

        if (p00 == p10 && p01 == p11) return 1;
        else if (p00 == p11 && p01 == p10) return -1;
      }
    }
    return (index)0;
}


/*
    The idea here is to match 3 points with the points
        of some face.  There are 3! == 6 ways (permutations)
        that the three point numbers can match the selected points.
*/

index points_are_a_face(ps0, ps1, ps2)
uindex ps0, ps1, ps2;
{
    uindex f, p0, p1, p2;

    for (f = 0; f < faces; f++)        /* search all faces */
    {
        p0 = face[f].p[0];
        p1 = face[f].p[1];
        p2 = face[f].p[2];

        if ((p0 == ps0 && p1 == ps1 && p2 == ps2) ||
            (p0 == ps0 && p1 == ps2 && p2 == ps1) ||
            (p0 == ps1 && p1 == ps0 && p2 == ps2) ||
            (p0 == ps1 && p1 == ps2 && p2 == ps0) ||
            (p0 == ps2 && p1 == ps0 && p2 == ps1) ||
            (p0 == ps2 && p1 == ps1 && p2 == ps0))
        {
            return (index)f;
        }
    }
    return (index)-1;
}


/*    Same as the above function, only it works with two points,
    returning a face that the two points (edge) share (if there is one).
    This function could be used to see if two points are "connected",
    or to find a face adjacent to an edge.  For this second usage,
    we don't want to find the original face, so we pass its number
    so that we can avoid it.
*/

index face_next_to_edge(ps0, ps1, avoid_f)
uindex ps0, ps1;
index avoid_f;
{
    uindex f, p0, p1, p2;

    for (f = 0; f < faces; f++)        /* search all faces */
    {
        p0 = face[f].p[0];
        p1 = face[f].p[1];
        p2 = face[f].p[2];

        if ((p0 == ps0 && p1 == ps1) ||
            (p1 == ps0 && p0 == ps1) ||
            (p0 == ps0 && p2 == ps1) ||
            (p2 == ps0 && p0 == ps1) ||
            (p1 == ps0 && p2 == ps1) ||
            (p2 == ps0 && p1 == ps1))
        {
            if (f != avoid_f) return (index)f;
        }
    }
    return (index)-1;
}


/* check if a point is a member of any face */

index point_is_on_a_face(ps0)
uindex ps0;
{
    uindex f;

    for (f = 0; f < faces; f++)        /* search all faces */
        if (face[f].p[0] == ps0 || face[f].p[1] == ps0 || face[f].p[2] == ps0)
            return (index)f;
    return (index)-1;
}


/* check if a point is a member of any poly */

index point_is_on_a_poly(ps0)
uindex ps0;
{
    uindex po, v;

    for (po = 0; po < polys; po++)        /* search all faces */
        for (v = 0; v < poly[po].verts; v++)
            if (poly[po].p[v] == ps0)
                return (index)po;
    return (index)-1;
}


index points_are_on_a_poly(array, verts)
uindex *array, verts;
{
    uindex po, pverts, test;
    register uindex i, j;

    if (!polys) return (index)-1;

    for (po = 0; po < polys; po++)
    {
        pverts = poly[po].verts;
        if (pverts < verts) continue;    /* skip simpler polys */

        for (i = 0; i < pverts; i++)
        {
            /* test forwards */
            test = 1;
            for (j = 0; j < verts; j++)
                if (array[j] != poly[po].p[(i+j) % pverts])
                    { test = 0; break; }
            if (test) return (index)po;

            /* test backwards */
            test = 1;
            for (j = 0; j < verts; j++)
                if (array[verts-j-1] != poly[po].p[(i+j) % pverts])
                    { test = 0; break; }
            if (test) return (index)po;
        }
    }
    return -1;
}


/*
    Like face_next_to_edge, except it works with polys.
    You can optionally pass a poly index for it to ignore.
    Returns -1 on failure.
*/

index poly_next_to_edge(ps0, ps1, avoid_p)
uindex ps0, ps1;
index avoid_p;
{
    uindex p, p0, p1, v, verts;

    for (p = 0; p < polys; p++)        /* search all polys */
    {
        verts = poly[p].verts;
        p0 = poly[p].p[verts-1];
        for (v = 0; v < verts; v++)
        {
            p1 = poly[p].p[v];

            if ((p0 == ps0 && p1 == ps1) ||
                (p1 == ps0 && p0 == ps1))
            {
                if (p != avoid_p) return (index)p;
            }
            p0 = p1;
        }
    }
    return (index)-1;
}


void zoom(buttons)
int buttons;
{
    register coord oldgx = gx, oldgy = gy, oldgz = gz;

    if (buttons & 1)
        SCALE /= 0.9;
     else
        SCALE *= 0.9;

    convert_to_xyz(mx, my);

    offset.x += oldgx - gx;
    offset.y += oldgy - gy;
    offset.z += oldgz - gz;
    convert_points();
    refresh_all();
}


void set_rot_axis(x, y, z)
coord x, y, z;
{
    axis.x = x;
    axis.y = y;
    axis.z = z;

    mode = M_ROTP;
    draw_bars();
    rot_angle = 0.0;
    draw_rot();    /* Redraw angle */
}


void rotate()
{
    register index i1;
    double  co = cos(rot_angle),
            si = sin(rot_angle);
    int move_p;

    if (flags.rotate_all)
        move_p = points;
    else
        move_p = selected;

    if (!enough_room(0, 0, 0, 0, 0, 0, move_p, 0))
    {
        /* this operation will overflow the undo buffer! */
        if (!tell_user_overflow()) return;
    }

    begin_operation();
    if (flags.rotate_all)
    {
        for (i1 = 0; i1 < points; i1++)
          rotate_point(i1, co, si);
    }
    else
    {
        for (i1 = 0; i1 < selected; i1++)
          rotate_point(select_p[i1], co, si);
    }
    end_operation();
    mode = M_ROTA;
    flags.copy_made = 0;
    refresh_all();
}


void rotate_point(p, co, si)
index p;    /* the point to rotate */
double co, si;    /* the cosine and sine of the angle to rotate */
{
    coord dx, dy, dz;

    dx = point[p].x - axis.x;
    dy = point[p].y - axis.y;
    dz = point[p].z - axis.z;

    switch (dmode)
    {
        case DM_XY:
            move_point_to(p, axis.x + (coord)(co*dx - si*dy),
                             axis.y + (coord)(si*dx + co*dy), point[p].z);
            break;
        case DM_XZ:
            move_point_to(p, axis.x + (coord)(co*dx - si*dz), point[p].y,
                             axis.z + (coord)(si*dx + co*dz));
            break;
        case DM_ZY:
            move_point_to(p, point[p].x, axis.y + (coord)(si*dz + co*dy),
                             axis.z + (coord)(co*dz - si*dy));
            break;
        case DM_3D:
            break;
    }
    convert_point(p);
}


void find_center(cx, cy, cz)
double *cx, *cy, *cz;
{
    long sum_x, sum_y, sum_z;
    register unsigned int i;

    *cx = *cy = *cz = 0;

    if (!points) return;

    sum_x = sum_y = sum_z = 0;
    for (i = 0; i < points; i++)
    {
        sum_x += point[i].x;
        sum_y += point[i].y;
        sum_z += point[i].z;
    }
    *cx = (double)sum_x / points;
    *cy = (double)sum_y / points;
    *cz = (double)sum_z / points;
}


void find_radius(cx, cy, cz, radius)
double cx, cy, cz, *radius;
{
    register unsigned int i;
    coord dx, dy, dz;
    double dist = 0.0;

    for (i = 0; i < points; i++)
    {
        dx = point[i].x - cx;
        dy = point[i].y - cy;
        dz = point[i].z - cz;
        dist += sqrt((double)(dx*dx + dy*dy + dz*dz));
    }
    *radius = dist / (double)points;
}


/*
    This integer sqrt function is about *15* times faster than
    the standard double-precision math function under lattice.
    It was the result of a thread on comp.graphics a while ago.
    I've tested the precise accuracy of it's results up to 600000.
     I suppose if i was a real CS god i could PROVE that this algo
    always works.
*/

unsigned long isqrt(v)
unsigned long v;
{
    register long t = 1L<<30, r = 0, s;

#define STEP(k) s = t + r; r >>= 1; if (s <= v) { v -= s; r |= t;}

    STEP(15);    t >>= 2;
    STEP(14);    t >>= 2;
    STEP(13);    t >>= 2;
    STEP(12);    t >>= 2;
    STEP(11);    t >>= 2;
    STEP(10);    t >>= 2;
    STEP(9);    t >>= 2;
    STEP(8);    t >>= 2;
    STEP(7);    t >>= 2;
    STEP(6);    t >>= 2;
    STEP(5);    t >>= 2;
    STEP(4);    t >>= 2;
    STEP(3);    t >>= 2;
    STEP(2);    t >>= 2;
    STEP(1);    t >>= 2;
    STEP(0);

    return (ULONG)r;
}


/* return fibonacci element (sum of natural numbers up to n) */

fibo(n)
{
    int sum = 0;

    while (n)
    {
        sum += n;
        n--;
    }
    return sum;
}


void set_clockwisdom()
{
    index f, q;

    for (f = 0; f < faces; f++)
    {
        if (determine_clockwisdom(face[f].p[0], face[f].p[1], face[f].p[2]))
        {
            change_face_clockwisdom(f);
        }
    }
    for (q = 0; q < polys; q++)
    {
        if (determine_clockwisdom(poly[q].p[0], poly[q].p[1], poly[q].p[2]))
        {
            change_face_clockwisdom(f);
        }
    }
}


int determine_clockwisdom(p0, p1, p2)
index p0, p1, p2;
{
    Point3D center;
    Vector3D v0, v1, v2, normal;
    coord dot;

    center.x = zero;
    center.y = zero;
    center.z = zero;

    v0.x = point[p0].x - center.x;
    v0.y = point[p0].y - center.y;
    v0.z = point[p0].z - center.z;

    v1.x = point[p1].x - center.x;
    v1.y = point[p1].y - center.y;
    v1.z = point[p1].z - center.z;

    v2.x = point[p2].x - center.x;
    v2.y = point[p2].y - center.y;
    v2.z = point[p2].z - center.z;

    /* now we find the normal by finding the cross product v0 X v1 */
    normal.x = v0.y * v1.z - v0.z * v1.y;
    normal.y = v0.z * v1.x - v0.x * v1.z;
    normal.z = v0.x * v1.y - v0.y * v1.x;
    /* note that there is no need to normalize the normal */

    /* we now have a definition of the plane through {center, p0, p1} */
    /* the dot product of v2 with the plane's normal gives the distance
        to the plane, of which we are only interested in the sign */
    dot = v2.x * normal.x + v2.y * normal.y + v2.z * normal.z;
    return (dot > zero);
}


void change_face_clockwisdom(f)
index f;
{
    index tmp;
    tmp = face[f].p[1];
    face[f].p[1] = face[f].p[2];
    face[f].p[2] = tmp;
}


void change_poly_clockwisdom(q)
index q;
{
    index tmp;
    index pi;    /* point index */
    int verts = poly[q].verts;

    for (pi = 1; pi < (verts+1)/2; pi++)
    {
        tmp = poly[q].p[pi];
        poly[q].p[pi] = poly[q].p[verts - pi];
        poly[q].p[verts - pi] = tmp;
    }
}